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ABSTRACT 

Context. Several efforts are currently made from space missions in order to get accurate solar astrometric measurements i.e. to probe 
the long term variations of solar radius or shape, their link with solar irradiance variations and their influence on earth climate. These 
space missions use full disk solar imagery. In order to test our ability to perform such measurements from ground on the long term, we 
need to use similar techniques and instruments simultaneously from ground and space. This should help us to model and understand 
how the atmosphere affect ground based metrologic measurements. However, using full imagery from ground instead of the traditional 
astrolabe technique immediatly raise the question of the effect of refraction and how well we can correct from it. 
Aims. The goal is to study in details the influence of pure astronomical refraction on solar metrologic measurements made from 
ground-based full disk imagery and to provide the tools for correcting the measurements and estimating the associated uncertainties. 
Methods. We use both analytical and numerical methods in order to confront commonly or historically used approximations and exact 
solutions. 

Results. We provide the exact formulae for correcting solar radius measurements at any heliographic angle and for any zenith distance. 
We show that these corrections can be applyed up to 80° of zenith distance provided that full numerical integration of the refraction 
integral is used. We also provide estimates of the absolute uncertainties associated with the differential refraction corrections and 
shows that approximate formulae can be used up to 80° of zenith distance for computing these uncertainties. For a given instrumental 
setup and the knowledge of the uncertainties associated with local weather records, this can be used to fix the maximum zenith 
distance one can observe depending on the required astrometric accuracy. 
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'- ■ 1 . Introduction order to use the same technique simultaneously from ground and 

^■"l space and in order to inter-calibrate the different measurements. 

^j- ; Ground based solar astrometric measurements have up to now Using full disk imager y f rom ground raise however the question 
O ■ been base d on t he so-called equal altitudes method (Debarbat of the influen ce of astronomical refraction and how well we can 
(f) & Guinot [1970l) . They have historically been made from tran- correct for it 
t— I ' sit instruments or astrolabes. Several instruments, derived from 

L" : Danjon astrolabe, have been dedicated to solar diameter mea- The effect of atmospheric refraction is to change the true 

surements, as DORaySol experiment (Morand et al. ESD- to P ocentnc zemth an 8 le z of a celestlal ob J ect t0 a lower ob " 
Observations consist in determining the transit times, through served one z - ^ refraction function R{z) is defined by: 
^ i the same equal zenith distance circle, of the two solar limbs 
which are the extremities of a vertical solar diameter. Accuracy 

of these measurements is then mainly determined by datation z — z' - R(z) (1) 

accuracy and not by the optical resolution of the instrument. As 
the two limbs are observed at equal zenith distances, influence of 

astron omical refraction is inherently reduced (e.g. Laclare et al. Alternatively, we may take the true angles as argument and de- 
II 9961 1 . Only the small climatic conditions variations (tempera- fine the associated refraction function R by: 
ture, pressure, relative humidity) between the two crossings, dis- 
tant from a few minutes of time, can still play a role. The main 

drawback of this method is that only vertical diameters can be z - z > _ R(^) (2) 

determined. 

Recent work in the field of solar metrology involve mea- 
surements from space using full disk solar images (Dame et al. If the refraction function R(z) is known, the associated function 
[19991 Kuhn et al. [20T2b . PICARD-SOL (Meftah et al. l20L3l) is R(z') can easily be evaluated for any true zenith distance z' by 
a ground based project that was set up at Calern observatory in solving the non linear equation x - R(z' - x) = 0. 
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2. fundamental equations for astronomical 
refraction 

From Snell's law of refraction applied to a spherical atmosphere, 
the curvature of a light path is linked to the local refractive index 
n through the so-called refractive invariant : 



n r sin(£) = constant 



(3) 



where ^ is the local zenith distance i.e. the angle between the 
light ray and the radius vector r from Earth center. From this, 
the differential refraction along the light ray is obtained by: 



dn 
dR = -tan£ — 

n 



(4) 



In order to find the total amount of refraction at observer posi- 
tion, we can integrate along the full ray path from n - « b s an d 
£, — z at observer position up to n — 1 outside the atmosphere. 

f* dn 

R = tan£ — (5) 

J\ n 

This can be done either by direct numerical integration of Eq . (0 
after an appropriate change of variable (Auer & Standish 2000) 
or by using a full ray-tracing procedure solving the system of 
coupled differential equations provided by Eq. (|4]i and the dif- 
ferentiation of Eq. (0 (van der Werf 2003 2008 1. This, in prin- 
ciple, requires a model of the full atmosphere i.e. temperature, 
pressure, density etc. .at any point through the light path. In the 
next section we recall why this is in fact not needed if we avoid 
areas close to the horizon and give some usual approximations 
of the refraction integral. 

3. Approximation to the refraction integral 

For zenith distance up to 70°, the refraction integral can be 
evaluated with good accuracy without any hypothesis about the 
structure of the atmosphere: it depends only on temperature and 
pressure at the observer (Oriani's theorem, see also: Ball 1908, 
Young 2004). This justifies that, over time, a large number of 
nearly equivalent approximate formulae have been derived that 
do not require the full knowledge of the structure of the real 
atmosphere. A development of the refraction integral into semi- 
convergent series of odd power of tan(z) is what is commonly 
found in textbooks (e. g. Ball [19081 Smart [19651 Woolard & 
Clemence 1966} Danjon [T980l ). An example of this will be given 
in Sect. 13.11 In fact the first two terms of such expansion (up 
to tan 3 ) corresponds to what is known as Laplace formulae of 
which Fletcher (U931| l said that no reasonable theory differs by 
more than a few thousandths, hundredths, tenths of a second at 
z = 60°, 70°, 75° respectively. 

For large zenith distance, tan(z) power series will diverge at 
the horizon and are not appropriate. Closed formula valid at low 
zenith distance and that are finite at the horizon can however 
still be found (see e.g. Wittmann 1997). Assuming an exponen- 
tial law for the variation of air density with height, it's possi- 
ble for instance to derive a formula involving the error function 
(Fletcher 1931 1 Danjon 1980). Another example is Cassini's ex- 
act formula for an homogeneous atmosphere model. While phys- 
ically un-realistic, the model of Cassini, thanks to Oriani's theo- 
rem, gives also excellent results up to at least 70° of zenith dis- 
tance while remaining finite down to the horizon (Young 2004). 
For large zenith distances however, Young (2004) have shown 
that the lowest layers of the atmosphere and especially the lapse 
rate at observer becomes progressively dominant as one observe 



closer to the horizon. This therefore should be included in atmo- 
spheric models and we can not avoid anymore the full numerical 
evaluation of the refraction integral. 

In the following sub-sections we present first in details the 
refraction model as it was used for reducing solar astrolabe data 
at Calern observatory, then we give the full error function model 
from which the Calern model was actually derived and finally 
we recall Cassini's formula. In Sect. [5] these three approxima- 
tions will then be compared to full numerical integration of the 
refraction integral using a standard atmosphere model. 

3. 1. Refraction model used at Calern observatory for Solar 
metrology 

The refraction model that was used for the reduction of astro- 
labe measurements at Calern observatory is a truncation of the 
expansion in odd power of tan(z) (Danjon 1980). For an observer 
at geodetic latitude <p and altitude h above the reference ellip- 
soid, the refraction R is obtained as a function of the observed 
zenith angle, the wavelength (A) and local atmospheric condi- 
tions i.e. pressure (P), absolute temperature (T), and relative hu- 
midity(/ A €[0,l])by: 



R(z, A, P, T,f h ,h,tp) = a( 1 - P) tan(z) - a(J3 - § ) tan 3 (z) 

\2 



+3a (/?-§) tan 5 (z) 



where 

a(T,P,f h ,A) = n obs - 1 



(6) 



(7) 



is the air refractivity for local atmospheric conditions and the 
given wavelength, and 



RT,h,<p) = €<T)/r l fy>,h) 



(8) 



is the ratio between the height I of the homogeneous atmosphere 
and the earth radius of curvature r c at observer position. The ho- 
mogeneous atmosphere has by definition a constant air density 
p equal to the one at observer position and its height is such that 
it would give the same pressure as the one recorded at observer 
position. Note that we do not assume here that the atmosphere 
is homogeneous, we just use the reduced height that can be ob- 
tained for any real atmosphere just from the pressure and density 
at observer. Assuming furthermore ideal gas law for dry air we 
have: 



t(T) 



P _ P o T 

Pg Po go T ' 



(9) 



where p = 1.293 kgm" 3 for T = 273.15 K, P = 101325 Pa 
and normal gravity go - 9.80665 ms~ 2 . The radius of curvature 
for Calern observatory ((p - 43°45'7", h = 1323 m) was approx- 
imated by the minimum reference ellipsoid curvature at latitude 
45° and sea level (Chollet 1981 see Appendix©: 



r e (45°,0) = 6367.512 km 



(10) 



Ambient air refractivity was deduced from the refractive index 
no(A) under standard conditions and the partial pressure of water 
vapor p by applying the formula recommended by the first reso- 
lution of the 13th General Assembly of the International Union 
of Geodesy and Geophysics (IUGG [19631 Baldini [19631 . After 
conversion to Pa (Pascal) as the pressure unit, the equation be- 
comes: 

a(T,P,f h ,A) = ^l(n (A)- 1)£-4.13 10~ 10 p(f h ,T)\ (11) 
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Refractivity under standard condition (sea level, T - To, P — Pq, 
0% humidity, 0.03% of carbon dioxide) was taken from the work 
of Barrel and Sears (IT9391 : 



n (A) -l = \ 2876.04 + 



16.288 0.136 



(10 6 T) 2 (10M) 4 



10" 



(12) 



Partial pressure of water vapor for the current temperature and 
relative humidity was deduced from a fit of water vapor pressure 
data published by the Bureau Des Longitudes (1975) for tem- 
peratures between -15 °C and +25 °C. The resulting equation, 
converted to Pa, is (Chollet fTWTT l: 



p(fh,T) 



■ f h 6.1075 10 2 e 7 - 2921(r2 ( 7 '- 7 '(>)-2.84i(r 4 (7--7-„) 2 



(13) 



Finally, we note that local atmospheric pressure P was measured 
from the height H (in mm) of a mercurial barometer and its 
temperature 9 (in °C). Taking into account corrections for lo- 
cal gravity (latitude and altitude) and for temperature (through 
the volume thermal expansion of mercury and the coefficient of 
linear thermal expansion of the tube), P was obtained b}(3 (see 
AppendixIBl: 

P = H[l -2.6410~ 3 cos(2^>)- 1.96 10 7 /i - 1.63 10 4 ^} (14) 



3.2. Error function formula 

In fact, in Eq. ||6}, only the first two terms which correspond to 
Laplace formula can be found without any hypothesis on the real 
atmosphere (only the reduced height I and the refractivity at ob- 
server are needed). The term in tan 5 comes from an additional 
assumption, namely the fact that air density follows an expo- 
nential decrease with height (actually with a well chosen vari- 
able which vary almost linearly with height, see Danjon (1980), 
Fletcher (1931)). This leads to the following equation: 



R = a 



2 - a 



sin(z) ¥ 



cos(z) 
V2y6 " a, 



(15) 



with : 
¥(x) = e 



"1°°^ 2 



e x (1 - erf(x)) 



(16) 



from which Eq. © was derived by keeping only the three first 
terms of its asymptotic expansion. 

3.3. Cassini 

By comparing the results with a full integration method, Young 
(2004) shows the superiority of Cassini's formula over the 
series-expansion approach and advocates its use by astronomers. 
Cassini assumed an homogeneous atmosphere for which he ob- 
tained the exact formula: 



. j n ohs r c sin(z)\ . /r c sin(z) 
R = asm - asm 



r c + e 



r c + £ 



(17) 



Again, it can be shown (Ball 1908) that expanding this formula 
also leads to the to first two terms of Eq. (O i.e. to Laplace for- 
mula. 




Fig. 1. Geometry for the solar shape due to astronomical refrac- 
tion. The dashed circle represents the true solar disk of centre C 
and radius Rq while the elliptical shape (full line) represents the 
observed Sun of centre C. The point at the top represents ob- 
server's zenith. 

4. On the observed shape of the Sun due to pure 
astronomical refraction 

In this section, we assume that the Sun is a perfect sphere of 
angular radius Rq at 1 AU and that there is no other effect af- 
fecting its observed shape than astronomical refraction defined 
by Eq. ©. 

In the horizontal coordinate system (zenith distance- 
azimuth), we note (z' , A©) the true position of the Sun cen- 
tre (C') observed at zenith angle z©; (z', A) the true position 
of a point (If) of the solar limb observed at zenith angle z; 
Sz - z - Zq and SA - A - A©. Figure Q] shows all the angles in- 
volved. Each true limb point position can be defined by the angle 
if/' e [-7r, 7r[ between the direction C'L' and the vertical circle. 
Similarly, each observed limb point can be located by the angle 
if/ G [-n, n[ between the observed direction CL and the vertical 
circle. However, because the figure is symmetric with respect to 
the vertical circle, we consider only the interval [0,7r] for if/ and 
if/' in the following. For observation with an Alt-Az mount this 
would correspond directly to the angle with one of the CCD axis. 
For an equatorial mount, one CCD axis is aligned with the hour 
circle passing through the celestial poles and the Sun and there- 
fore the vertical circle can be materialized on the solar image by 
computing first the parallactic angle between these two circles. 

If d(ifr) = d{if/') is the angular distance between the observed 
position of the Sun centre and the observed limb points, we de- 
fine by: 



<d>- 



i r 

n Jo 



d(if/) 2 dif/ 



1/2 



i r 

n Jo 



d(if/fdif/ 



1/2 



(18) 



the geometric mean radius of the observed Sun. The horizontal 
and vertical angular extent of the observed Sun are noted D/, and 
D v respectively and, the flattening is given by: 



/ = 



D h 



(19) 



Chollet ( 1981 ) used erroneously 2.64 10 4 in this equation. 
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Following Mignard (2002), we define the magnification F as the 
ratio between the vertical size of the image (Sz) of a small object 
to its true size (&')■ From Eqs. (Q} and (0, we have: 



= dz _ ! _ dR _ (df\ = / rfR\" 
~ rfz' ~ <fe { ~ 1 dz I ~ \ + dz J 



(20) 



The distorsion A is then defined as the rate of change of the mag- 
nification: 



A- — -_r 2 — - -- — 
~Tz~~ ~dz T ~~Td7i 



4.1. Approximate formulae for all zenith angles 



(21) 



Any limb point true position can be located by its projections on 
the vertical circle passing through the true Sun centre, and on the 
great circle perpendicular to this vertical circle passing through 
the limb point (see Fig. [TJ. Because all the angles involved are 
small, we can write: 



x' = RQCOs(ift') 
y< = R Q sin(^) 

and: 



(22) 
(23) 



x' + y = R 



O 



By looking at the expression of the observed values x and y of 
these projections, one can obtain an approximate formula for the 
observed shape of the Sun. 

The projection x' on the vertical circle can be approximated 
by keeping the two first terms of a Taylor expansion of the re- 
fraction: 



x* » z f - z'q = Sz + R(z) - R(zq) * Sz (l + y ) 



^> 



The observed projection y is linked to z and 6 A both by the cosine 
and sine rules: 



cos(y') =* cos 2 (z f ) + sin 2 (z') cos(c)A) 
sin(y) = sin(&4) sin(z') 



(26) 
(27) 



Differentiating Eq. (l26l l and using Eq. (l27l i with sin(y') » y', 
sin(SA) m SA and dz' - -R(z') leads to: 



dy' = 



-y r(z') 



tan(z') 
The observed distance y is then obtained by: 

R(z') 



y ^ SA sin(z) = / + dy' = y ( 1 - 



tan(z') 



(28) 



(29) 



Finally, by reporting Eqs. ((25) and (|29l in Eq. (l24l > and using 
Eqs. (l20b and (f2TT >. we obtain: 



Sz _ A /SzY 

r 2 IrJ 



5A sin z 



1 - 



tan(z') 



= R, 







(30) 



where the magnification and distortion are taken at zq. From 
this we can deduce the position of the two vertical limb points 



and the observed vertical extent of the image. For A <K Rq and 
SA - 0, we find: 



d(n) =* T^o 1 + 
d(0) * TRq\1 



ARq 

2 
ARq 



(31) 
(32) 



and thus: 

D v =d(0) + d(n)^2TR Q (33) 

In the horizontal direction we obtain from Eq. (l30l with Sz = 0: 



D h = 2d(7r/2) a 2Rq 



1 - 



tanfeJJ 



(34) 



4.2. Approximate formulae for small zenith angles - elliptic 
shape 

Keeping only the first term in Eq. (0 is equivalent to neglecting 
Earth curvature. We obtain the following approximation valid 
close to the zenith only {z < 45°): 



R(z) = fctan(z) with : k = a(l-0) 



(24 s ) ^ or tn ^ s fl at -Earth approximation we can also write: 



(35) 



(36) 



R(z') * k' tan(z') with : k' = k{\ - k sec 2 (z')) 

In that case and if we neglect the distortion, Eq. (l30t is reduced 
to the equation of a simple ellipse (see also e.g. Ball 1908): 



(l-£'sec 2 (z^)) 2 d"^) 2 



= R 







(37) 



where x - Sz and y - sin(z) SA can be assimilated to Cartesian 
coordinates on two perpendicular axes on the image. The major 
axis of the observed ellipse is thus given by: 



y=*0(l-*') 

while the observed minor axis is: 



^=tf (l-*'sec 2 (z^)) 



(38) 



(39) 



We note from these equations that the Sun is shrunken in all di- 
rections. The observed horizontal diameter is smaller than the 
true diameter but remains the same for all zenith angles (c.f. 
Fig|7]i while the observed vertical diameter decreases with in- 
creasing zenith distance. The combination of these two effects 
leads to the apparent flattening of the setting Sun (but keeping 
in mind that this approximate formula is not valid close to the 
horizon). From Eqs. ( [T9| >, (l38l and ( |39l , the flattening for small 
zenith angles is: 



/ - &tan 2 (4). 



(40) 



while, near the horizon, Eq. (l30l implies that the flattening is 
simply given by the vertical magnification taken at the the Sun's 
centre. For small zenith angles, the observed elliptic shape can 
be written as: 



d(ffr) 



D v 



: A /l-(2/-/ 2 )sin 2 (.A) 



(41) 
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which can be approximated by: 

d(if/) *Rq(\-K(i + cos 2 (if/) tan 2 (4))) , 

and the mean radius is obtained by: 



< d >= -I-J!_* » R Q 1 - k' - - tan 2 ^) 



(42) 



(43) 



more easily obtained using the sine rule rather than Eq. (146) . 
With sin(d(>r/2)) =* rf(7r/2), sin(fl©) =* /?© and cos(Rq) - 1, we 
obtain: 



d(>r/2) - sin(z) sin(M) = sin (z*, - £ (z^)) -^- 



4.3. Exacf formulae for all zenith angles 

The classical approximate formulae above are useful for under- 
standing the shape of the observed Sun in terms of magnifica- 
tion and distortion induced by refraction. Equation (l30l shows 
that the general shape is a distorted ellipse with more flattening 
in the lower part than in the upper's. However, the shape of the 
observed Sun can also easily be obtained, in the general case, 
without any approximation. In the following, we obtain first the 
solution of the forward problem: for given true Sun radius Rq 
and true zenith distance zL, we obtain the shape of the observed 
Sun for any given refraction model. Then, we give the solution of 
the inverse problem: from the observed solar shape, the knowl- 
edge of zk (from ephemeris) and assuming a refraction model, 
we deduce the true angular solar radius. 



4.4. Forward problem 

Here we assume that the true zenith distance of the Sun centre 
z'p, and its true angular radius Rq are known. For any refraction 

model R(z!), and true angle if/', we deduce the observed angle if/ 
and angular distance d(if/). Applying the cosine and sine formu- 
lae respectively, we have : 



z = acos 



6A = asin 



[cos (zq) cos (Rq) + sin (z^) sin (Rq) cos(^')] 

sin(J? Q ) sinfrfrV 

sin(z') 



(44) 



From Eq. (0, we can get the observed zenith distances: 
z = z'-R(z') and z Q = z' Q - R(z Q ) 



(45) 



and finally angular distances d(if/') between the observed Sun 
centre and the observed positions of each limb point are obtained 
by application of the cosine rule: 

d(ifr') -d(if/)- acosf cos(z) cos(z0) + sin(z) sin(z Q ) cos(&4) ) (46) 

where the observed angle if/ can be deduced from the true one by 
applying the sine rule: 

. /sin(&4)sin(z)\ . / sin(z) sin(fl ) . \ 
if/ = asm — ,-„.,^ | = asin I ■_,_,-, ,-„..^ sm w ) I (47) 



sinW)) 



sin(z') s in(d(^0) 



The smallest observed diameter of the Sun is obtained on the 
vertical direction: 

D v = d(0) + d(n) = 2fl - (R (z'q +Rq)-R (z' q - Rq)) (48) 

and the largest angular extent, observed in the direction parallel 
to the astronomical horizon is obtained by: 



D h = 2d(n/2) 



(49) 



We note that Eqs. ( l44b and (|43T > lead back to the approximation 
Eq. (l34l for the largest observed angular extent. This is however 



(50) 



■©' 



which, with a first order expansion of the sine function around 
Zq, leads to Eq. ( T34b . 

4.5. Inverse problem 

Here we give the solution of the inverse problem: given a re- 
fraction model (R(z), R(z')), knowing zL from ephemeris and 
the observed angular distance d(ifr) between the observed Sun 
centre and a limb point at an observed angle iff with the vertical 
circle, we deduce the true angular radius Rq. One can compute 
successively: 



"O 



SA 



<■© 



atan 



R(z'q) 



sin(i/>) 



sin(zQ)cot(<f (i/0) + cos(zq) cos(i/0 



z = asin 



sin(^) sm(d(if/)) 
sin(&4) 



z + R(z) 



(51) 



if/' - atan 



Rq = asin 



sin(&4) 

cos(zq) cos(&4) - sin(zj 

sin(&4) sin(z') 
sin(ifj') 



z©)cot(z f ) 



5. Results 

5.1. On the absolute value of refraction 

We first look at the absolute value of refraction and compare the 
various approximate formulae of Sect. |3]to the full numerical in- 
tegration of the refraction integral using a standard atmosphere 
(Sinclair 1982). This atmosphere is assumed to be spherically 
symmetric, in hydrostatic equilibrium and made of a mixture of 
dry air and water vapor that follows the perfect gas law. It is 
made of two layers: the troposphere with a constant temperature 
gradient which extends from the ground up to the tropopause 
at 1 1 km, and an upper isothermal stratosphere. Like in the US 
Standard Atmosphere (1976 1, the temperature and pressure at 
the surface are 288.15 K and 101325 Pa and the constant tropo- 
spheric lapse rate is 6.5 K km -1 . In the troposphere, the relative 
humidity //, is assumed constant and equal to its value at the ob- 
server. The partial pressure of water vapor in a tropospheri c layer 
at temperature T is then obtained by: 



Pifh, T) 



fh (247.1) 



10 2 



(52) 



which, with 6 = 18.36, never depart by more than 0.5 hPa from 
Eq. ( fT3l for temperature lower than 30°. Dry air is assumed in 
the stratosphere. Finally, Eq. ([TO and its derivatives with respect 
to T and P are used to find air refractivity along the integral path. 
The numerical integration was performed by using the 
method of Auer & Sandish (2000) also recommended by the 
Astronomical Almanac (Seidelmann 1992). The program used 
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zenith distance [deg] 



Fig. 2. Absolute differences (in mas) between a reference model 
and the different approximate refraction formulae as a function 
of the true zenith distance. The reference model is obtained 
by full numerical integration of the US Standard Atmosphere 
( 1 19761 1 and Ciddor (119961 1 equation for air refractivity. From 
top to bottom: tan 5 expansion Eqs. d6t-([T3li. full error function 
Eq. ( fT3l l. tan 5 expansion Eq. ©, Cassini's formula Eq. (fTTT i. All 
approximate formulae but the top one use Ciddor (11996) air re- 
fractivity 



is based on the one published by Hohenkerk & Sinclair (1985 1 
but adapted in order to use a dispersion equation based on the 
work of Peck & Reeder (11972) in replacement of the less accu- 
rate equation of Barrel & Sears (IT939) (Eq. (H2)). For T = 15 °C, 
P - Po, 0% humidity and 0.045% of carbon dioxide, we take: 



n (A) - 1 = 



0.05792105 



0.00167917 



238.0185 -(10 6 ^K 57.362 - (10<U)~ 



(53) 



This dispersion equation was also used by Ciddor (1996 1 who 
derived a new set of equations for calculating the refractive in- 
dex of air which was subsequently adopted by the International 
Association of Geodesy (I AG 1999) as a new standard. In the 
following, all computations have been made using A = 535.7 nm 
which is one of the wavelengths used by the PICARD-SOL 
project. 

Figure [2] shows the absolute differences in milliarcseconds 
(mas) between the approximate formulae and the exact integral 
evaluation for zenith distances up to 80°. We immediately see 
that for zenith distance lower than 75°, all the approximate for- 
mulae lead to less than 50 mas of absolute error. The full line 
corresponds to the tan 5 formula Eq. (0 described in Sect. 13.11 
while the dashed line corresponds to the same formula but using 
the new Ciddor (1996) equations instead of Eqs. dTTT >-(fT3l> for 
computing air refractivity. For zenith distances lower than 80°, 
the impact of using the old formula for refractivity never exceed 
80 mas. The superiority of Ciddor equations to better fit observa- 
tions and this for a wider range in wavelengths is however clearly 
established. The two other lines correspond to the error function 
(dot-dash) and Cassini (triple dots-dash) formulae both using the 
Ciddor ( 1996) equation for refractivity. These two last formulae 
were selected mainly because, unlike the series expansions in 
tan(z), they are finite at the horizon. The full integration with 
standard atmosphere conditions leads to a refraction of about 
1980" at the horizon. The error function and Casini formulae 
lead respectively to 2088" and 1180" corresponding to relative 
errors of 5% and 40% respectively. This tends to favour the use 
of the error function formula over Cassini's one very close to the 
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Fig. 3. Relative error on refraction as a function of zenith 
distance for different tropospheric lapse rate. The reference 
model use US Standard Atmosphere (1976| with a lapse rate 
of 6.8 KknT 1 . The top curve correspond to an isothermal model 
and other atmosphere models have lapse rate of 2.5, 5, 10 and 
7.5 KkrrT 1 (from top to bottom at low zenith distance). All mod- 
els are computed using full numerical integration. 



horizon. The hypothesis made to derive the error function for- 
mula are indeed more realistic than Cassini's hypothesis of an 
homogeneous atmosphere. It has however been shown that re- 
fraction below 5° of the horizon is variable and strongly depend 
on the local lapse rate and properties of the boundary layer above 
or below the observer's eye (e.g. Young 2004l. Within few de- 
grees from the horizon, refraction may be influenced by thermal 
inversion boundary layers, ducting or other phenomena leading 
to extreme refraction. In this range, the local lapse rate must be 
known and it is not expected that any formula using just the tem- 
perature and pressure at observer could give an accurate absolute 
refraction. 

It is however probably more interesting to look in the range 
between 60° and 85° of zenith distance, which is more impor- 
tant to astronomers willing to push in that range the limits of 
their astrometric measurements using only temperature and pres- 
sure recorded at observer position. We first note from Fig. 
that, between 60° and 77°, the tan 5 expansion formula is ac- 
tually giving slightly better absolute refraction values than the 
error function formula. If we now assume that temperature and 
pressure at observer position are perfectly known, the only re- 
maining important unknown in the atmospheric model is the tro- 
pospheric lapse rate. We can however fix limits for a realistic 
lapse rate: it must lie between an isothermal model and a lapse 
rate of 10 KknT 1 which would correspond to an adiabatic atmo- 
sphere (Young 2004). Figure [3] shows the absolute value of the 
relative error for such models with lapse rate ranging from to 
10 K knT when they are compared to the standard model with a 
lapse rate of 6.8 K knT 1 . From this we can deduce that, no mat- 
ter what is the real atmosphere, if the conditions at observer are 
known, the relative error on refraction is lower than 0.01% for 
zenith angles below 77° and lower than 0.4% for zenith angles 
between 77° and 85°. 

5.2. On the mean solar radius correction 

Figure |4] shows the difference between the true radius of the Sun 
and the mean radius of the observed Sun as defined by Eq. (fT8l 
as a function of the true zenith distance of the centre of the Sun. 
The exact formula Eq. ( f46b was used and we took standard con- 
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Zenith distance [deg] 



V [deg] 



Fig. 4. Difference between the true solar radius and the 
observed one as a function of the true zenith distance. 
The full line corresponds to average weather conditions at 
Calern (T = 15 °C, P = 875 hPa). The dashed and dot-dashed 
lines correspond respectively to T — -10 °C, P - 900 hPa and 
T — 30 °C, P = 850 hPa. All calculations are made using the ex- 
act formulae Eqs. ( fT8l and d46i i for Calern station assuming 50% 
humidity. 



Fig. 6. Difference between the true solar radius and the angu- 
lar distances between the observed Sun centre and the observed 
positions of each limb points between the vertical (north for 
if/ = 0° and south for if/ = 180°) and the horizon (if/ = 90°). The 
full lines are for z' Q = 70°, 50°, 30° and 10° respectively from 
top to bottom and are for average weather conditions at Calern. 
The dashed and dot-dashed lines are for zk = 70° and the same 
extreme weather conditions as in Fig. |4] 




Fig. 5. Difference between the correction due to refraction on 
the mean solar radius as calculated from integrating the exact 
formula Eq. ( l46l l or using the approximate formula Eq. (l43l . The 
dashed line is obained by replacing k' by k in Eq. (l43l 




true zenith distance [deg] 



Fig. 7. Contraction of the horizontal radius (Rq - d(n/2)) as a 
function of the true zenith distance zL. The full line is for av- 
erage weather conditions at Calern. The dashed and dot-dashed 
lines are for the same extreme weather conditions as in Fig. |4] 



ditions for Calern observatory (T = 15 °C, P = 875 hPa). The 
dashed and dot-dashed lines are for T = -10 °C, P = 900 hPa 
and T — 30 °C, P - 850 hPa respectively in order to illustrate 
the maximum amplitude of the effect at Calern station. The 
difference in the mean radius correction between the two ex- 
treme weather conditions range from 50 mas at the zenith up 
to 1850 mas at z' = 85°. It reaches 100 mas around z' — 55° and 
200 mas around z' - 70°. This represents always less than 0.2% 
of the correction. 

Figure [5] shows the difference between the exact formula 
obtained by integrating Eq. (|46] > and the approximate formula 
Eq. d43l corresponding to an elliptical shape. The dashed line 
illustrates the result if k' is approximated by k (see Eq. d36i>). In 
both cases the difference remains less than 20 mas for zenith dis- 
tances lower than 70°. For larger zenith distances however, errors 
increase rapidly and the refraction function should be evaluated 
using full numerical integration. 



5.3. On the angular dependence of solar radius correction 

For precise metrologic measurements of the Sun and in order 
to correct for other effects (optical aberrations, turbulence, etc..) 
that are dependent on the position on the image, one may want 
to correct not the mean radius but each individual radius mea- 
sured at all angles if/. This can be done by following the proce- 
dure given in Sect. 14.51 Figure [6] is obtained from Eq. (l46t and 
illustrates the amplitude of the correction as a function of if/ for 
different values of z%., the true zenith distance of the Sun cen- 
tre. We see that the horizontal diameter (if/ - 90°) is affected by 
refraction (by about 2 x 0.23" = 0.46" for the chosen weather 
conditions) in agreement with Eq. ( l34l . The north and south ver- 
tical corrections (if/ = 0° and 180° respectively) are also slightly 
different in agreement with Eqs. (l3T1l-(l32l. Figure Q shows that 
the contraction of the horizontal radius lies between 210 and 
260 mas depending on the actual weather conditions and re- 
mains constant for all zenith distances below 80°. It then de- 
creases rapidly towards zero at the horizon. 
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true zenith distance [cleg] 



Fig. 8. Partial derivatives of the vertical diameter correction (£,,) 
as a function of the true zenith distance. Partial derivatives in 
temperature, pressure, zenith distance and relative humidity are 
given in mas KT 1 , mas hPa 1 , mas arcmin -1 and mas/% from top 
to bottom (at 40°) respectively. The full, dashed and dot-dashed 
lines are for the same weather conditions as on Fig. |4] 




true zenith distance [deg] 



Fig. 9. Uncertainties on the vertical diameter correction as- 
suming AT - 0.5 K (dotted line), AP - 1 hPa (dashed line), 
A//, = 5% and Az' Q = 5.4', 1.0' or 0.2' (full lines from top to 
bottom). The total error is obtained by summing the four contri- 
butions. 

5.4. On uncertainties associated to radius corrections 

We have shown that, apart from weather conditions at observer's 
position, differences in atmospheric models and especially dif- 
ferent tropospheric lapses rate will not play any significant role 
at least up to 85° of zenith distance. The four main contribu- 
tions are therefore uncertainties in temperature, pressure, humid- 
ity and, for large zenith distance, uncertainties on the true zenith 
distance itself. 



dd(tf/) 



dXi 



AdW^^AX, X = {T,P,f h , Z ' Q } 



(54) 



It should be noted that we assume here observations made using 
filters with a narrow bandwidth around A. For broadband filters, 
an additional term dd(tf/)/dA should be added by differentiating 
Eq. 1531 . The largest uncertainty will be obtained for the ver- 
tical diameter (D v - d(Q) + d(ii)) which is the most affected by 
refraction. Figure [8] shows the four partial derivatives contribut- 
ing to AD V between the two extreme weather conditions chosen 
above for Calern (see Sect. 15. 2\ . The partial derivatives shown 



have been obtained by numerically differentiating Eq. d46b but 
we have also checked that the analytical expressions that can be 
derived from the approximate elliptical shape Eq. (l42t are actu- 
ally valid up to 80° of zenith distance. Closer to the horizon the 
partial derivative over the zenith distance becomes significantly 
overestimated (c.f. Fig.fTOt. From Eq. (l42l and taking k' ^ k, we 
obtain: 



C dd(tfr) 
dXi 


= 


dk 
dX t 


ddffl 

I «4 




2k co 



(l+cos 2 (if,)t<m 2 (z' Q ))Ro j=1..3 



= 2k cos 2 (i/r)sec 2 (sQ tan (z' Q ) Rq 



and from Eqs. ©, ©, (TTT), ([35}. (g3, we obtain: 



' & = -C l ^(no(A)-l) + C 3 T s - i f h (c 2 - 6 -^-) 



(55) 



dk 

BP 



en _ f (t-\ 



C 2 )(n Q (A) - 1) 



(56) 



M- = -C 3 (T- 1 - C 2 )T 5 
v ofh 



where: 



d = T /P , C 2 l = drcpogo, C 3 = 4.13 10- 8 r () (247.ir> (57) 

For temperature, pressure and humidity, we assume uncertain- 
ties of AT = 0.5 K, AP = 1 hPa and Af h = 5% which are typical 
for a standard weather station. The precision on the true zenith 
distance relies on ephemeris calculations and a correct timing. 
At any given time ephemeris can give not only zL but also the 
instantaneous rate dzU,/dt and, from the knowledge of the image 
exposure time At, one can deduce an uncertainty on zL by: 



Az' 



■0 



d-' 







dt 



At 



(58) 



The maximum rate is about 650" min -1 at summer solstice. 
Image exposures of 1 s, 5.5 s or 30 s would then correspond 
to a maximum uncertainty Azk of 0.18', 1' or 5.4' respectively. 
Figure|9]shows the contribution of these uncertainties to the total 
uncertainty on vertical diameter correction for large zenith dis- 
tances. We can see for instance that for 1 ' precision on the zenith 
distance (or 5.5 s exposure), the uncertainty coming from zenith 
distance can become, above 70°, of the same importance as the 
combined uncertainties coming from temperature and pressure 
records. The total relative error on the vertical diameter cor- 
rection (AD v /(2Rq - D v )) remains however below 1% up to 



*© 



85°. 



6. Conclusions 

We have obtained in Sec. 14.51 the exact formulae that can be used 
to correct solar radius measurements at any heliographic angle 
and any zenith distance from the effect of astronomical refrac- 
tion. Using full integration of the refraction integral, we have 
shown that this correction can be applied for any true zenith dis- 
tance up to 85° with a relative uncertainty on the correction that 
remains less than 1%. Absolute uncertainties on these correc- 
tions are also derived that allows us to fix the maximum zenith 
distance one should observe depending on the needed metro- 
logic accuracy. Figure [TO] shows the maximum total absolute 
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uncertainty obtained on the solar radius assuming that the ver- 
tical radii have been observed at different zenith distances. We 
use A = 535.7 nm and two exposure times used by the PICARD- 
SOL project at this wavelength. Because we took the maximum 
value for dzL/dt, this curves represent only upper limits, the ac- 
tual value of dzL/dt should be use for each measurement. From 
this, one can deduce that observing below 70°, 75° or 80° of 
zenith distances will keep the absolute uncertainties on refrac- 
tion corrections below 10, 20 and 50 mas respectively. The com- 
parison between numerical derivatives (full lines) and the use of 
approximate formulae Eqs. (I54t-(l58l (dashed lines) shows that, 
even if the approximate formulae should not be used above 70° 
for correcting the measurements (c.f. Fig. |5}, they can be used at 
least up to zk = 80° for estimating the uncertainties. 

In summary, the process that we suggest to correct ground 
based radii measurements from refraction for true zenith dis- 
tances up to 80° is as follow. Inputs are: the measurements d(<p) 
and eventually their associated errors Sd(<p) where (f> is an arbi- 
trary angle defined on the solar image; the time of image record 
and the exposure time At; weather records (P, T, fh) and their 
associated uncertainties (AT, AP and A//,); the wavelength (A) 
and observer's geodetic coordinates (<p, h). One can then succes- 
sively: 

- find the direction of the zenith on the image and associate 
each angle <p to its corresponding angle ijj (c.f. Fig. []]). 
Depending on the instrumental setup, this may require the 
computation of the parallactic angle from ephemeris, 

- determine z' and dzL/dt from ephemeris at the time of im- 
age record, 

- calculate Rq using Eqs. ( BIT ) and full numerical integration 
for the refraction function R(z, A, P, T, //,, h, ip), 

- estimate Ad(t//) from Eqs. (I54ll-(l5"8ll and the knowledge of 
AT, AP, Af h , At and dz' Q /dt, 

- estimate ARq from: 



AR =R< 



Ad(ip) + 6d(tf/) 
~dW) ' 



0- 



(59) 



For zenith distances lower than 70° full numerical integration 
can be replaced by Eq. © in order to evaluate the refraction 
function (c.f. Fig. [2}. In both cases Ciddor (1996) equations 
should be used for computing air refractivity at observer po- 
sition. The corresponding codes are available from the authors 
upon request. 

It is important to keep in mind that, at all zenith distances, 
other phenomena such as extinction or optical turbulence must 
be taken into account for ground based solar metrology. We 
know that they will dominate refraction effects at low zenith dis- 
tances. Close to the horizon extinction is proportional to refrac- 
tion (Laplace's extinction theorem) and effects of optical turbu- 
lence (e.g. Ikhlef et al. 2012 and reference therein) will become 
increasingly important knowing that the Fried parameter varies 
as sec(z) -06 . It is interesting however to know that for any zenith 
distance up to 80° refraction can be reliably corrected and un- 
certainties on this correction estimated. After these correction 
are applied, all other phenomena impacting metrologic measure- 
ments can therefore be investigated without fearing contamina- 
tion by astronomical refraction even at high zenith distances. 
The mean radius correction presented here (c.f. Fig. as well 
as mean turbulence corrections have been applied to correct the 
first PICARD-SOL measurements (Meftah et al. 120131) . The cor- 
rections that can be applied individually for each heliographic 




true zenith distance [deg] 



Fig. 10. Maximum absolute uncertainties on radius esti- 
mates using measurements in the vertical direction corrected 
from refraction. This curves are obtained for /l = 535.7nm, 
T = (15 ± 0.5) °C, P = (875 ± 1) hPa, /,,=50%±5% and 
dz'^/dt = 650" min -1 . The top and bottom lines are for 8.9 s 
and 1.3 s of exposure time respectively. Full lines give the 
results from full numerical derivatives calculations while dashed 
lines are obtained using approximate formulae Eqs. (I54t- (l58b . 

angles should be used in future work in order to disentangle the 
different effects. 

Finally we note that we have considered only the radial 
symmetric-component of refraction also called pure or normal 
refraction. There also exists an asymmetric component known as 
anomalous refraction (e.g. Teleki 1979 1 resulting from the tilted 
atmospheric layers. Anomalous refraction may depend not only 
on zenith distance but also on azimuth and it can lead to seasonal 
or high frequency effects (see e.g. Hirt (2006) and references 
therein). The amplitude of such effect has however been found to 
be lower than 0.2" for local effects and one order of magnitude 
less for regional effects that may originate higher in the atmo- 
sphere (e.g. Hu 1991). Moreover it has been shown that anoma- 
lous refraction is spatially coherent at scales of at least 2° (Pier 
et al. 2003) and it has been established from dedicated observa- 
tions that its main source is confined in the layer immediately 
above ground level (less than 60 m, see Taylor et al. (2013)). It 
is therefore difficult to believe that differential effects of anoma- 
lous refraction and especially the one that may be triggered in 
the Upper Troposphere - Lower Stratosphere (UTLS) interface 
(c.f. Badache-Damiani et al. 2007) could lead to significant bias 
on solar astrometric measurements relying on direct solar disk 
imaging. 
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Appendix A: Note on the radius of curvature at 
Calern observatory 

According to the WGS84 reference ellipsoid, the Earth's equa- 
torial and polar radii are given respectively by a = 6378 . 137 km 
and b = 6356.752 km. The curvature in the (north-south) merid- 
ian and at the geodetic latitude of Calern solar astrometric in- 
struments tf = 43°45'7" is then given by: 

(ab)2 — (A.1) 



r° = 



(a 2 cos 2 (ip) + b 2 sm 2 (ip)) 



3/2 



= 6365.985 km 



One could also consider the mean radius of curvature calculated 
for Calern. From the curvature in the prime vertical (normal to 
the meridian): 



r?° = 



J a 2 cos 2 (</j) + b 2 sin 2 (i^) 



= 6388.371 km 



(A.2) 



we can deduce the radius of curvature for any azimuth angle A 
by: 

— : ?— (A.3) 

cos 2 (A) sin"(A) 



^ 



from which we can deduce the mean radius of curvature averag- 
ing over all directions, by: 



<r r > 



= Jfrf = 



a-b 



a 2 cos 2 (<^>) + b 2 sin (ip) 



= 6377.168 km (A.4) 



If, instead of the radius of curvature, one considers the distance 
from geocenter, we have: 



R 



la 4 cos 2 (<p) + b 4 sin 2 (<^) 
a 2 cos 2 ((p) + b 2 sin 2 (<^) 



6367.955 km 



(A.5) 



One should add to these values the elevation of the observer 
above the reference ellipsoid (h = 1.323 km for Calern observa- 
tory). If we consider that, on average, we observe the sun closer 
to the north-south direction than east- west direction we can take: 



r c = r" + /i = 6367.308 km 



(A.6) 



which is very close to the value used by Chollet (119811 1. 

Finally we note that, for ephemeris calculations, the 
geodedic latitude should be corrected for the local gravimetric 
deflection. For Calem solar astrometric instruments this lead to 
an astronomic latitude tp ast = 43°44'53" which is also compat- 
ible within 1" with the direct measurements made using a full 
entry pupil astrolabe on the same site. Similarly, we note that 
taking into account the local undulation with respect to the refer- 
ence ellipsoid leads to a height above sea level of h s \ = 1.271 km 
for Calern solar astrometric station. 

Appendix B: Note on the corrections applied to 
mercurial barometer reading 

The two corrections (for gravity and barometer temperature) can 
be written as multiplicative factors (e.g. Princo 12007) : 



P = H 



1 +L( 

1 + M 8) go 



(B.l) 



where P is the corrected atmospheric pressure, H is the barom- 
eter reading, M = 1.818 10~ 4 K is the coefficient of volume 
thermal expansion of mercury, and L = 1 .84 10~ 5 K" 1 is the co- 
efficient of linear thermal expansion of brass. According to the 
1967 reference system formula (Helmert's equation), we have: 



g - g45 (l - acos(2ip) - bcos 2 (2tfi)) 



(B.2) 



where #45 = 9.8061999 ms is the gravity acceleration at mid 
latitude, a = 2.64 10~ 3 and b = 1.96 10~ 6 . This can be corrected 
from the so-called Free Air Correction (FAC) which accounts 
for the fact that gravity decreases with height above sea level 
(Cfac = -3.086 10~ 6 s~ 2 ), itself corrected in order to take into 
account the increasing gravity due to the extra mass assumed for 
a flat terrain (Bouger correction, Cb = 4.2 10~ 10 m 3 s~ 2 kg -1 ). 
For a mean rock density of p r = 2.67 10 3 kg trT 3 this leads to: 



C g = (C FAC +P,-C b ) = -1.96 10- 6 s- 



(B.3) 



Close to 45° of latitude, the second term of Eq. (IB. 2b can be 
neglected and, if we note e = 1 — g4s/go = 4.6 10 s , Eq. (IB. II ) 
can be approximated by: 

P=H(1- e) (l - (M-L) 6) jl -acos(2<p) + — h\ (B.4) 

Neglecting second order terms leads to Eq. (fT4l . 
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We note that absolute gravity measurements have now 
been made at Calern geodetic observatory leading to 
g = (980215549.2 ± 12.6) lO^ms" 2 (Nicolas et al. 120061 
This shows that the relative error on the correction g/go dis- 
cussed above and previously used for refraction calculations 
was less than 5 10~ 5 . One could however now directly use 
Eq. ( IB. II ) with the measured value of local gravity. 
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